#include<stdio.h>
#include<stdlib.h>
#include<math.h>

/* Retorna o valor máximo absoluto de um array de floats */
float Max(float array[], int tam){
	float max = fabs(array[0]);
	int i;

	for(i = 0; i < tam; i++)
		if(fabs(array[i]) > max)
			max = fabs(array[i]);

	return max;
}

int main(){
	int i, j, k, J_ORDER, J_ROW_TEST, J_ITE_MAX;
	float J_ERROR, diag, soma = 0;
	float resA = 0, resB, iter_final, erro;
	float **MA, *VB;
	float **MA2, *VB2;
	float **x, *diff;

	/* Leitura dos dados */
	scanf("%d", &J_ORDER);
	scanf("%d", &J_ROW_TEST);
	scanf("%f", &J_ERROR);
	scanf("%d", &J_ITE_MAX);

	/* Alocação de memória das variaveis */
	MA = (float**)malloc(J_ORDER * sizeof(float*));
	MA2 = (float**)malloc(J_ORDER * sizeof(float*));
	for(i = 0; i < J_ORDER; i++){
		MA[i] = (float*)malloc(J_ORDER * sizeof(float));
		MA2[i] = (float*)malloc(J_ORDER * sizeof(float));
	}

	x = (float**)malloc((J_ITE_MAX+2) * sizeof(float*));
	for (i = 0; i < J_ITE_MAX+2; i++) 
		x[i] = (float*)malloc(J_ORDER * sizeof(float));

	VB = (float*)malloc(J_ORDER * sizeof(float));
	VB2 = (float*)malloc(J_ORDER * sizeof(float));
	diff = (float*)malloc(J_ORDER * sizeof(float));

	/* Leitura da matriz MA */
	for(i = 0; i < J_ORDER; i++)
		for(j = 0; j < J_ORDER; j++)
			scanf("%f", &MA[i][j]);

	/* Leitura do vetor VB */
	for(i = 0; i < J_ORDER; i++)
		scanf("%f", &VB[i]);

	/* Criacao de MA* e VB* (neste caso MA2 e VB2) */
	for(i = 0; i < J_ORDER; i++){
		diag = MA[i][i];
		MA2[i][i] = 0;
		for(j = 0; j < J_ORDER; j++){
			if(j != i)
				MA2[i][j] = MA[i][j]/diag;
			else
				MA2[i][j] = 0;
		}
		VB2[i] = VB[i]/diag;
	}

	/* Calculo da matriz x e do erro */
	for(i = 0; i < J_ORDER; i++)
		x[0][i] = VB2[i];
	for(k = 0; k < J_ITE_MAX; k++){
		for(i = 0; i < J_ORDER; i++){
			soma = 0;
			for(j = 0; j < J_ORDER; j++)
				soma += MA2[i][j] * x[k][j];
			x[k+1][i] = VB2[i] - soma;
			diff[i] = fabs(x[k+1][i] - x[k][i]);
		}
		erro = (Max(diff,J_ORDER) / Max(x[k+1],J_ORDER));

		if(erro <= J_ERROR)
			break;
	}

	if(k == J_ITE_MAX)
		k--;

	/* Calculo da resposta */
	for(i = 0; i < J_ORDER; i++)
		resA += MA[J_ROW_TEST][i] * x[k+1][i];
	resB = VB[J_ROW_TEST];

	printf("Iterations: %d\n", k+1);
	printf("RowTest: %d => [%f] =? %f\n", J_ROW_TEST, resA, resB);

	/* Libera memória */
	free(VB);
	free(VB2);
	free(diff);

	for (i = 0; i < J_ORDER; i++) {
  		free(MA[i]);
		free(MA2[i]);
	}
	free(MA);
	free(MA2);

	for (i = 0; i < J_ITE_MAX+2; i++)
		free(x[i]);
	free(x);

	return 0;
}
